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The phase diagram of star polymer solutions in a good solvent is obtained over a wide range of 
densities and arm numbers by Monte Carlo simulations. The effective interaction between the stars 
is modeled by an ultrasoft pair potential which is logarithmic in the core-core distance. Among the 
stable phases are a fluid as well as body-centered cubic, face-centered cubic, body-centered orthog- 
onal, and diamond crystals. In a limited range of arm numbers, reentrant melting and reentrant 
freezing transitions occur for increasing density. 

PACS numbers: 64.70.-p, 82.70.Dd, 61.25.Hq 
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A major challenge in statistical physics is to under- 
stand and predict the macroscopic phase behavior from 
a microscopic many-body theory for a given interaction 
between the particles 0] . For a simple classical fluid Q , 
this interaction is specified in terms of a radially sym- 
metric pair potential V(r) where r is the particle sepa- 
ration. Significant progress has been made during the 
last decades in predicting the thermodynamically stable 
phases for simple intermolecular pair potentials, such as 
for Lennard-Jones systems, plasmas or hard spheres, us- 
ing computer simulations 0J and density functional the- 
ory |jj . An important realization of classical many-body 
systems are suspensions of colloidal particles dispersed 
in a fluid medium. A striking advantage of such col- 
loidal samples over molecular ones is that their effective 
pair interaction is eminently tunable through experimen- 
tal control of particle and solvent properties 0). This 
brings about more extreme pair interactions, leading to 
novel phase transformations. For instance, if the colloidal 
particles are sterically stabilized against coagulation, the 
'softness' of the interparticlc repulsion is governed by the 
length of the polymer chains grafted onto the colloidal 
surface, their surface grafting density and solvent qual- 
ity. Computer simulations and theory have revealed that 
a fluid freezes into a body-centered-cubic (bec) crystal 
for soft long-ranged repulsions and into a face-centered- 
cubic (fee) one for strong short-ranged repulsions ||. 
This was confirmed in experiments on sterically stabi- 
lized colloidal particles j(| . A similar behavior occurs for 
charge-stabilized suspensions where the softness of V(r) 
is now controlled by the concentration of added salt . 
Less common effects were observed for potentials involv- 
ing an attractive part aside from a repulsive core. In 
reducing the range of the attraction, a vanishing liquid 
phase has been observed || and an isostructural solid- 
solid transition was predicted H . More complicated pair 
potentials can even lead to stable quasicrystalline phases 



and a quadruple point in the phase diagram 1 1 1 



The aim of this letter is to study the phase diagram 
of an ultrasoft repulsive pair potential V(r) which is lo- 
garithmic in r inside a core of diameter a and vanishes 
exponentially in r outside the core. The motivation to 
do this is twofold: first, such a potential is a good model 
for the effective interaction between star polymers in a 
good solvent jll],[l2]], which can be regarded as steri- 
cally stabilized particles where the size of the particles 
is much smaller than the length of the grafted polymer 
chains JL3|| . These stars are characterized by their arm 
number (or functionality) /, i.e., the number of polymer 
chains tethered to the central particle, and their corona 
diameter a which measures the spatial extension of the 
monomer density around a single star center. Second, 
more fundamentally, phase transitions for such soft po- 
tentials are expected to be rather different from that for 
stronger repulsions. From a study of the pure logarith- 
mic potential in two spatial dimensions jl4j , it is known 
that one needs a critical prefactor to freeze the system, 
which is quite different from, e.g., inverse-power poten- 
tials. Furthermore, the potential crossover at r — a is 
expected to influence drastically the freezing transition, 
if the number density p of the stars is near the overlap 
concentration, p* w l/er 3 . 

We obtain the full phase diagram of star polymer solu- 
tions by Monte Carlo simulation and theory. As a result, 
among the stable phases are a fluid as well as bee, fee, 
body-centered-orthogonal (bco), and diamond crystals. 
We emphasize that the stability of a bco crystal with 
anisotropic rectangular elementary cell and a diamond 
structure was never obtained before for a radially sym- 
metric pair potential. In fact, there is a widespread belief 
in the literature that anisotropic or three-body forces are 
solely responsible for a stable diamond lattice [^5|. We 
show that both the crossover at r — a and the ultrasoft- 
ness of the core are crucial for the stability of the bco and 
the diamond phase. Moreover, we get reentrant melting 
for 34^/^60, and reentrant freezing for 44^/^60 as p 
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is increasing. Some features of the presented phase dia- 
gram have already been observed in a syatem of copoly- 
mer micelles exhibiting a very similar interaction to star 
polymers |6[[il|. 

With ksT denoting the thermal energy, our effective 
pair potential between two star centers is a combination 
of a logarithm inside the core of size a and a Yukawa- 
potential outside the core p^ ]: 



V(r) = ^k B Tf> 



-Hi) 



cxp(-y7(r- ff )/2 t r) 



< 0; 

(r > a), 
(1) 



such that both the potential and its first derivative (or, 
equivalently, the force) are continuous at r = a . The 
decay length of the exponential is given by the largest 
blob diameter within the Daoud-Cotton theory for single 
star polymers Jl7[ . Experimental support for this poten- 
tial comes from neutron scattering data on the structural 
ordering of 18-arm stars in the fluid phase [^2) and shear 
moduli measurements in the crystalline phase of micelles 
JTsfl . Furthermore, microscopic simulations of two star 
polymers have shown that this potential provides an ex- 
cellent description of the effective star interaction for a 
broad range of arm numbers jlfit] . We note that V(r) 
becomes the hard sphere potential for / — > oo. 

Due to the purely entropic origin of the interstar repul- 
sion, the strength of the pair potential ([]]) scales linearly 
with ksT, causing the temperature to be an irrelevant 
thermodynamic quantity. Therefore, for the calculation 
of the phase diagram, only the packing fraction of the 
stars, 77 = 7r/6/9cr 3 , and the arm number / matter, the 
latter playing the role of an "effective inverse tempera- 
ture." We use computer simulations to access the phase 
diagram. The free energies of the fluid phase and several 
possible solid phases are calculated by thermodynamic 
integration via Monte Carlo simulations pc[ . The free 
energy of the fluid phase, Ffi, is obtained either by the 
well-known "pressure- or density- route" or, alter- 

natively, by the so-called "/-route" . The pressure-route 
relates the free energy for nonvanishing 77 to that at zero 
packing fraction, keeping / fixed. In the /-route, / is 
used as an artificial thermodynamic variable, now keep- 
ing 77 fixed. The free energy of star polymers with a 
certain arm number / is then obtained by the following 
integration: 
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Here, U = Yli<j ^(l r i — r jl) ^ s the total potential energy 



function which depends on / since V(r) depends on / 
parametrically. (...) /' denotes the canonical ensemble av- 
erage for a system with fixed arm number /'. Therefore, 
in order to carry out the /-route integration, a series of 



simulations at fixed 77 but for increasing /' is performed 
to calculate the integrand of Eq. (||) . 

We use the Frenkel-Ladd method for continuous poten- 
tials to obtain the free energy of the solid phases |^0|,|l| . 
For these Monte Carlo calculations, suitable candidate 
crystal structures have to be chosen. Our method to 
get information about the possible stable structures for 
fixed / and 77 consists of two steps: first, we calculate 
lattice sums for a wide class of crystals, including the 
'usual' structures with cubic elementary cells (fee, bec, 
hep, and simple cubic) and several 'unusual' structures. 
These unusual structures are the hexagonal lattice, the 
diamond lattice, representations of quasicrystalline struc- 
tures (see, e.g., Ref. @), and generalizations of the usual 
structures, which were obtained by stretching the ele- 
mentary cell lengths (denoted as a, b, and c) of these 
structures by arbitrary factors, then using the two in- 
dependent ratios b/a and c/a as minimization parame- 
ters of the lattice sum. Second, we calculate the global 
bond order parameters |p2f of the equilibrated structures, 
which were spontaneously formed in a first set of simu- 
lations, always starting from a purely random configu- 
ration. Crystal structures whose bond order parameters 
are in agreement with these measured parameters, and 
which have reasonably small values of the lattice sum, are 
then chosen as candidate structures for the free energy 
calculations. This procedure was performed for a wide 
range of arm numbers, 18 < / < 512, and packing frac- 
tions < 77 < 1.5. Finally, the obtained free energy data 
at fixed / were used to explore the phase boundaries via 
the common double tangent construction. The resulting 
phase diagram is displayed in Fig. [I]. 
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FIG. 1. The phase diagram of star polymer solutions for 
different arm numbers / versus packing fraction 77. The 
squares and the circles indicate the phase boundaries as ob- 
tained from computer simulations and theory, respectively; 
lines are only guides to the eye. The statistical error of the 
simulations is of the order of the symbol size. The triangles 
indicate the freezing and melting point of hard spheres. 

In the explored range of / and 77, four different stable 
crystal structures are found besides a fluid phase. For 
/ < f c ~ 34, the fluid phase is stable for all densities, 



2 



Typeset using REVTeX 



which is in agreement with results obtained from a effec- 
tive hard sphere mapping procedure p3j ] and from scaling 
theory [[fl|. We remark that Witten et al. (ll]] only es- 
timated f c within one order of magnitude to be around 
/ ~ fOO. For / > f c , at least one stable crystal phase is 
found. We focus first on the crystal phases at 0.2^77^0.7: 
for f c < /^54, a bcc phase is found, whereas for /^70, 
only the fee structure turns out to be stable. At interme- 
diate / (54^/^70), bcc- fee phase transitions occur. For 
0.2^7^0.7, the mean interparticle distance r = p -1 / 3 is 
larger than a, leaving only the exponential part of V{r) 
to be relevant for the phase behavior. Therefore, the 
observance of a fee phase for large /, corresponding to a 
short-range, strongly screened potential, and a bcc phase 
for small /, corresponding to a long-range, less screened 
potential, is analogous to the phase behavior found for 
charged colloids d,Q] . In Fig. [l], the freezing and melting 
points for hard spheres, corresponding to / — * 00, are 
shown as well, denoted by black triangles. We empha- 
size that even star polymers with very high arm num- 
bers freeze at considerably smaller 77 than hard spheres. 
In fact, our simulations show that a 'hard-sphere like' 
structure is only found for extremely high arm numbers 
10000. Thus the change in the phase boundary can- 
not be shown on the scale of the figure. 

Let us now consider the phase behavior for 77 ~ 0.7, 
where F is in the order of a and the logarithmic part of 
V(r) becomes relevant. From our calculations, a reen- 
trant melting transition, i.e., a transition from a solid to 
a liquid phase with increasing 77, is found for 34 < /^60. 
We note that this reentrant melting was already pre- 
dicted qualitatively by Witten et al. pj. For /^60, a 
solid-solid phase transformation into a bco phase takes 
place. This unusual phase is stable up to 77 ~ 1.0. For 
44^/^60, the remolten liquid refreezes into this bco 
structure at 77 ~ 0.80. At 77 ~ 1.0, a further solid-solid 
phase transition from the bco into a diamond structure 
is found, the latter being stable for arm numbers / ^44 
and packing fractions up to 77 ~ 1.4-1.5. Notice that the 
extension of the two phase regions ("density-jumps") of 
all encountered phase transitions is extremely small due 
to the soft character of V(r) p3L Moreover, the empiri- 
cal Hansen- Verlet freezing rule |25[| is valid for all points 
at the phase boundaries where we calculated the static 
structure factor S(q). This also includes the reentrant 
melting transition for 77 ~ 0.7, where the S(q) for the 
fluid begins to show unusual behavior [p3| . 

We develop now a physical intuition for the unusual 
occurrence of the bco and diamond phase. For this pur- 
pose, we report first on the detailed structure of the bco 
phase. At fixed 77, the bco crystal is described by the two 
lengths ratios of its elementary cell, b/a and c/a, respec- 
tively. In order to calculate the free energy of the bco 
crystal by the Frenkel-Ladd method, these ratios had to 
be determined from a first set of simulations. In these 



jVpT-simulations [|20|, the system was free to adopt its 
optimal values for b/a and c/a, starting either from a 
purely random configuration or a initial bco configura- 
tion. Within the error bars, the so determined elemen- 
tary cell length ratios were in agreement with the values 
obtained from the minimization of the lattice sums. We 
therefore took the lattice sum results as input for the 
free energy calculations. These ratios increase with 77 
from b/a ~ 2.24 and c/a ~ 1.32 at 7/ = 0.7 to b/a ~ 3.14 
and c/a ~ 1.81 at 77 = 1.0 and are nearly independent of 
/. Fig. illustrates the resulting structure. 




FIG. 2. (a) A snapshot of a typical bco configuration for 
/ = 64 and r\ = 0.8 in a periodically repeated (cubic) sim- 
ulation box. The diameter of the spheres is the corona di- 
ameter <t; (b) Same as (a), now seen from the 'left' side of 
the simulation box shown in (a). Notice the high anisotropy 
of the lattice spacings. The elementary cell length ratios are 
b/a ~ 2.70 and c/a ~ 1.57. 

As can be seen from this figure, the anisotropy of the 
elementary cell leads to a strong interpenetration of the 
particle coronas along one of the lattice axes. In fact, 
over the whole stability range of the bco phase, the next 
neighbor distance along this axis is considerably smaller 
than a, whereas all other next neighbor distances are 
larger than a. This can be intuitively understood from 
the form of the potential (0) : due to the weak divergence 
for small r, there is no huge energy penalty in bringing 
the nearest neighbors close together. On the other hand, 
the potential falls off rapidly for r > a, so all the remain- 
ing neighbor shells are not costly in energy, too. With 
increasing 77, the distance of the two nearest neighbors in 
the bco is decreasing until the energy penalty becomes 
significant. Hence, the bco will then lose against another 
structure with more than two nearest neighbors inside the 
corona. A suitable structure is the diamond phase which 
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possesses four tetrahedrally ordered nearest neighbors. 
Indeed, our simulations show that all the other neighbors 
are kept outside the corona in the stability range range of 
the diamond. Therefore, both the ultrasoft logarithmic 
part and the crossover at r = a are crucial for the stabil- 
ity of the bco and the diamond phase. This provides a 
simple reason why such phases have not been found ear- 
lier for strongly repulsive interactions. We further note 
that the presented scenario also nicely expresses itself in 
the angle-average radial distribution functions g(r) of the 
bco and diamond solid, which show a similar anomaly as 
found in the g(r) of the fluid phase |25| ]. 

As for a further theoretical investigation, we solved the 
accurate Rogers- Young closure 1 26 to obtain the free en- 
ergy of the fluid for / = 18,32,40,48 and 64. For the 
aforementioned solid structures, we used the Einstein- 
crystal perturbation theory J2?J to calculate the associ- 
ated free energies. As this theory provides only an upper 
bound to the free energy, the domain of stability of the 
fluid is enhanced in comparison to the simulation results. 
The theory predicts 40 < f c < 48 and eliminates the do- 
main of stability of the bec crystal. Otherwise, as also 
shown in Fig. |l|, the same phase behavior as determined 
from simulations emerges. 

We finally note that all our predictions for p < 2p*, 
i.e. 77 ^ 1.0, should be verifiable in scattering experi- 
ments, since for these densities pair interactions are dom- 
inant. In fact, in recent experimental work on spherical 
diblock copolymer micelles, Cast and coworkers have al- 
ready confirmed a part of our results |^,^6| ■ The freezing 
transition in fee and bec crystals depending on the num- 
ber of arms / is found Q as well as reentrant melting 
with increasing rj jl(|. For the "most starlike" system, 
also a reentrant freezing is observed as predicted in Fig. 
[jj. For ry^ 1.0 however, when three stars exhibit over- 
laps within their coronae, many body interactions be- 
come important, which we have neglected in our calcu- 
lations using the pair potential (0). Nevertheless, from a 
theoretical point of view, this potential turned out to be 
interesting also for 77^ 1.0, resulting for the first time in a 
stable diamond structure for a purely radially symmetric 
pair interaction. 

In conclusion, we have determined the phase diagram 
of star polymers over a broad range of arm numbers / and 
packing fractions rj by computer simulations and theory. 
The phase diagram includes a fluid phase as well as four 
stable crystal phases. These crystal phases are a fee crys- 
tal and a bee crystal, as well as an unusual anisotropic 
bco structure and a diamond crystal. 
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